// @HEADER
//*********************************************************************//
//  SiYuan: A numerical PDE solver                                     //
//  Copyright (2022) YUAN Xi                                           //
//  This Software is released under the BSD 2-Clause license detailed  //
//  in the file "LICENSE" in the top-level SiYuan directory            //
//*********************************************************************//
// @HEADER

#ifndef _BC_NEUMANN_HPP
#define _BC_NEUMANN_HPP

#include "Teuchos_ParameterList.hpp"

#include "Panzer_STK_Interface.hpp"
#include "Panzer_GlobalIndexer.hpp"
#include "Panzer_PhysicsBlock.hpp"
#include "PanzerCore_config.hpp"
#include "Panzer_EvaluatorsRegistrar.hpp"
#include "Panzer_LinearObjContainer.hpp"
#include "Panzer_WorksetDescriptor.hpp"

#include "valfunctor.hpp"
#include "NeumannDef.hpp"

#include <string>
#include <map>

namespace SiYuan {

//**********************************************************************
template<typename EvalT, typename Traits>
class NeumannEvalutor : public PHX::EvaluatorWithBaseImpl<Traits>
{
  using ScalarT = typename EvalT::ScalarT;

  public:
    NeumannEvalutor( const NeumannBoundary& bc, const Teuchos::RCP<const panzer::GlobalIndexer> & indexer
      , Teuchos::RCP<panzer::PhysicsBlock>& side_pb );

  //  void postRegistrationSetup(typename panzer::Traits::SetupData d, PHX::FieldManager<panzer::Traits>& fm);
    void preEvaluate(typename Traits::PreEvalData d);
    void evaluateFields(typename Traits::EvalData d);

  protected:
    //! Finds the basis for the corresponding dof_name in the physics block.
    Teuchos::RCP<panzer::PureBasis> getBasis(const std::string dof_name,
					     const panzer::PhysicsBlock& side_pb) const;

    //! Allocates and returns the integration rule associated with an integration order and side physics block.
    Teuchos::RCP<panzer::IntegrationRule> buildIntegrationRule(const int integration_order,
							       const panzer::PhysicsBlock& side_pb) const;
                     
    Teuchos::RCP<panzer::LinearObjContainer> m_GhostedContainer;
    Teuchos::RCP<const panzer::GlobalIndexer> globalIndexer_;

  private:
    std::string m_residual_name;
    std::string m_flux_name;

    std::string  m_sideset_name;
    std::string  m_eblock_name;
    std::string  m_dof_name;
    std::string  m_bc_type;
    XYZLib::variable<ScalarT>  m_value;

    bool _need_scatter;

    // intermidiate result
    Kokkos::DynRankView<ScalarT, PHX::Device> neumann;
};


//**********************************************************************
template<typename EvalT, typename Traits>
class CLoadEvalutor: public PHX::EvaluatorWithBaseImpl<Traits>
{
//  using ScalarT = typename EvalT::ScalarT;

  public:
    CLoadEvalutor( const std::vector< CLoad >& bcs, const Teuchos::RCP<const panzer::GlobalIndexer> & indexer);

  //  void postRegistrationSetup(typename panzer::Traits::SetupData d, PHX::FieldManager<panzer::Traits>& fm);
    void preEvaluate(typename Traits::PreEvalData d);
    void evaluateFields(typename Traits::EvalData d);

  protected:

    //! Allocates and returns the integration rule associated with an integration order and side physics block.
    Teuchos::RCP<panzer::IntegrationRule> buildIntegrationRule(const int integration_order,
							       const panzer::PhysicsBlock& side_pb) const;
                     
    Teuchos::RCP<panzer::LinearObjContainer> m_GhostedContainer;
    Teuchos::RCP<const panzer::GlobalIndexer> globalIndexer_;

  private:
    std::string m_residual_name;
    std::string m_flux_name;

  //  std::string  m_nodeset_name;
    std::string  m_dof_name;
    std::string  m_bc_type;
    std::vector<panzer::LocalOrdinal> m_ldofs;
    std::map< panzer::LocalOrdinal, double > m_cloads;
    int called;
    //XYZLib::variable<ScalarT>  m_value;
};

}

#include "Neumann_impl.hpp"

#endif
